Fluid Flow Complexity in Fracture Networks: Analysis with Graph Theory and LBM

Through this research, embedded synthetic fracture networks in rock masses are studied. To analysis the fluid flow complexity in fracture networks with respect to the variation of connectivity patterns, two different approaches are employed, namely, the Lattice Boltzmann method and graph theory. The Lattice Boltzmann method is used to show the sensitivity of the permeability and fluid velocity distribution to synthetic fracture networks’ connectivity patterns. Furthermore, the fracture networks are mapped into the graphs, and the characteristics of these graphs are compared to the main spatial fracture networks. Among different characteristics of networks, we distinguish the modularity of networks and sub-graphs distributions. We map the flow regimes into the proper regions of the network’s modularity space. Also, for each type of fluid regime, corresponding motifs shapes are scaled. Implemented power law distributions of fracture length in spatial fracture networks yielded the same node’s degree distribution in transformed networks. Two general spatial networks are considered: random networks and networks with “hubness” properties mimicking a spatial damage zone (both with power law distribution of fracture length). In the first case, the fractures are embedded in uniformly distributed fracture sets; the second case covers spatial fracture zones. We prove numerically that the abnormal change (transition) in permeability is controlled by the hub growth rate. Also, comparing LBM results with the characteristic mean length of transformed networks’ links shows a reverse relationship between the aforementioned parameters. In addition, the abnormalities in advection through nodes are presented. characteristics of flow propagation like permeability [69]. The first steps in analyzing the fractures’ structural configurations included the fractality of the fracture systems. Having a fractal dimension with a small variation indicates that the system has a universality property. However, there is not such universal dimensionality in fracture networks. Researchers have recognized two fractal dimensions for the roughness (non-directional) of a single fracture [10]. During the past decade, developments in graph theory [11-15] have opened a new chapter in the study of large interwoven systems (complex systems). This new perspective tackles another facet of the complexities embedded in fractured materials, called structural or topological complexity. In contrast, statistical methods investigate uncertainties such as probability or, in more advanced forms, other uncertainties approaches (like fuzzy set theory, rough set theory or evidence theories). The relevant topics in this field might be summarized as follows: 1) transformation of spatial networks in graph forms; 2) information propagation (fluid flow, energy or waves) through networks; 3) deformation of fractures due to mechanical or hydro-mechanical forces and, more generally, the reaction of a disordered system (like a fracture system) to external forces. Our focus in this study is to generate different configurations of spatial fracture networks. In this pursuit, we try to use the simplest algorithms to fracture network generation. Then, we transform spatial fracture networks into graph forms where we ignore the spatial distribution of fractures. The advantage of these transformations is that it links the regular fracture networks to modern graph theory. The fracture zones as a high density of fractures also as the abnormal emergence of fractures are another part of our research. We present a simple algorithm that gives the directionality and effective radius of the hubs, while taking into account the growth of the hubs under a background growth of random joints. Advection of information (here fluid flow) through the generated spatial network and transformed networks is the main part of this study. To analyze fluid flow, three different methods were employed. We use lattice Boltzmann (LBM) and finite element methods (FEM) to obtain the fluid flow patterns, velocity, pressure distribution and permeability for the spatial fracture networks. Additionally, an advectionbased network equation is used for transformed networks. The organization of the paper is as follows. The first section presents a simple algorithm to generate fracture networks and show how they are converted into the graphs. In this section, we also introduce some basic characteristics of graphs. The next section summarizes three methods— LBM, FEM and advection-based networks— to model fluid flow (laminar) in fracture networks. The main part of our study will be covered in section 4, where we present the results and discuss the accuracy of the employed methods. Finally, the summary and conclusion of the present work is presented. 2. FRACTURE NETWORKS AND GRAPH THEORY In this part, a simple algorithm is introduced to cover the main characteristics of fracture networks. The transformation of constructed fracture networks into graph forms, and the distinguishing characteristics of the graphs, are demonstrated next. 2.1. Random Fracture Networks Several algorithms have been proposed to generate and cover the main statistical properties of natural fracture networks (based on the second approach presented in our introduction). The simplest algorithms in 2D consider the distribution of fracture length, dip direction of joint sets and joint spacing [7-9]. New generations of fracture networks algorithms (or Discrete Fracture Networks: DFN) consider other parameters that increase the accuracy of the estimations, such as: the fracture density parameter (number of fractures per length/area or volume) and fractal dimension, where the fracture center’s fractality is imposed through a hierarchal multiplicative process [9]. Most of the aforementioned algorithms (and codes) impose a power law or modified power law distribution to fracture length. The cut-off value in power law distribution shows a collapsing data set, indicating a similar mechanism/signature in the modeled data set. The basic idea behind the power law distribution of fractures remains a challenging question. However, recent developments in graph theory have tried to provide a reasonable answer. In three-dimensional scenarios, each individual fracture is assumed as triangle, circle or any other polygon shapes where the consecutive generations based on the pre-set parameters are accomplished [17]. Another main attribute of the fracture system is the density or spatial density of fractures. This property yields the possible concentration of fractures in space as well as fracture tips, around tunnels, the interface between two layers and generally any sharp or disturbed sub-fields within the system. The equivalent definition of such fracture zones may be found in modularity and the “hubness” concept, in which groups, communities or dense clusters, through random links, communicate with each other [18]. Our algorithm captures the power law distribution, directionality of the joint sets and the “hubness” properties of the networks, as described below (Figure 1). This algorithm is based on power law distribution of fracture length ( ) p l l γ − ∼ where l is the fracture length and gamma is the power that controls the connectivity and length. When gamma increases, the probability of finding fractures with long length reduces. Other parameters are: the number of joint sets and fracture zones (hubs), hub growth, background rock joint growth (hereafter called the external links or back joints), maximum diameter of growth in x and/or y directions, azimuths of joint sets and mean aperture of each fracture. The parameter related to “hub” and “back” growth is reported as reciprocal values. Then large values of growth show slow propagation of links. The spatial distributions of hubs are accomplished by employing a Gaussian distribution. The growth in hub zones and complementary background fractures are based on a power law distribution of fracture length and the uniform distribution of joint sets’ dips. We show that, for a wide range of gamma parameters, the distribution of links—in the transformed shape of fracture networks— obeys the power law distribution, indicating scale free networks. In other words, with the described method, we achieved modular networks with power distribution in transformed networks. n_g=200; %Number of Generations n=300; %Number of grids NJS =2; %Number of Joint set gama=.55; %parameter related to broadness of fracture length alpha=5*n./2; %parameter to scale the power law hub_growth=3;% reciprocal of growth rate of fractures around hubs Back_growth=2;%reciprocal of growth rate of background fractures ( ) ; 0.55 p l l γ γ − = ∼ N.F.Z=1:3 %number of fracture zones (hubs) r1 = round(10 + 10.*rand(1,1)); r2 = round(10 + 10.*rand(1,1)); %radius of hubs -effective zone of fracture zones (homogeneous cores) Fig. 1. Part of the parameters employed in fracture network code-The code was compiled in M-file format/Matlab. 2.2. Modern Graph Theory A network consists of the nodes and edges connecting them [11-12]. In order to analyze the network properties of the generated fracture networks, each fracture generation was mapped into a corresponding node in a graph space. When two fractures are intersected, two nodes are connected by a link. This procedure is illustrated in Figure 2, which enables us to investigate the networks using the tools from modern network theory discussed in the following. Fig. 2. The transformation procedure of a spatial fracture network in a graph shape with nodes and edges (also see [16]). Let us introduce some properties of the networks, such as: the clustering coefficient (C), degree distribution ( ( ) P k ) and average path length ( L ). The clustering coefficient describes the degree to which k neighbors of a particular node are connected to each other. Neighbors are defined as the connected nodes to a particular node. The clustering coefficient shows the collaboration between the connected nodes, i.e., the local structures. Assuming that the th i node has i k neighboring nodes, there can exist at most ( 1) / 2 i i k k − edges between the neighbors. i c is defined as the ratio: